"""
文件名称: ccp_user_association.py
功能: 使用CCP凸凹规划求解用户关联策略(内层)
作者: 喻越
创建日期: 2025-01-15
wtb修改与 2025-03-19
"""
from ..parameter.constant import const
import cvxpy as cp
import numpy as np
from .performance_calculations import compute_CE
from .system_initialization import *
from .performance_calculations import bs_energy_consumption_calculation


def function1(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req, r_bs_ue):
    ue_num = res_x_bs_ue.shape[0]
    c = const.P_WORK * bs_status + const.P_SLEEP * (const.ONE_BS - bs_status)  # shape =[1,bs_num]
    RB_W = np.tile(const.RB_W, (ue_num, 1))
    res = np.sum(RB_W * rb_req * r_bs_ue * res_x_bs_ue) * 3600 - ce_temp * const.ETA_GRID * np.sum(
        np.maximum((const.RB_PR / const.ETA) * np.expand_dims(np.sum(rb_req * res_x_bs_ue, 0), 0) + c - p_available,
                   const.ZERO_BS))
    return res


def function2(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req):
    c = const.P_WORK * bs_status + const.P_SLEEP * (const.ONE_BS - bs_status)  # shape =[1,bs_num]
    res = ce_temp * const.ETA_PV * np.sum(
        np.minimum((const.RB_PR / const.ETA) * np.expand_dims(np.sum(rb_req * res_x_bs_ue, 0), 0) + c, p_available))
    return res


def function(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req, r_bs_ue):
    return function1(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req, r_bs_ue) - function2(res_x_bs_ue, ce_temp,
                                                                                                bs_status, p_available,
                                                                                                rb_req)


def concave_derivation(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req):
    bs_num = res_x_bs_ue.shape[1]
    ue_num = res_x_bs_ue.shape[0]
    derivation = np.zeros(shape=[ue_num, bs_num])
    p_bs, _ = bs_energy_consumption_calculation(bs_status, res_x_bs_ue, rb_req)  # shape=[1,bs_num]
    for j in range(bs_num):
        for i in range(ue_num):
            if p_bs[0][j] < p_available[0][j]:
                derivation[i][j] = ce_temp * const.ETA_PV * (const.RB_PR / const.ETA) * rb_req[i][j]
            else:
                derivation[i][j] = 0
    return derivation


def optimize_concave(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req, r_bs_ue):
    bs_num = res_x_bs_ue.shape[1]
    ue_num = res_x_bs_ue.shape[0]
    bs_status_temp = np.tile(bs_status, (ue_num, 1))  # 复制拓展矩阵
    one_ue = np.ones(shape=[ue_num, 1])
    vec = np.array([0, 1])

    x_bs_ue = cp.Variable((ue_num, bs_num))
    # TODO 需要在这里添加约束，使得每个用户必有所属基站
    # cp.sum(x_bs_ue, 1, keepdims=True) == one_ue,

    constraints = [cp.constraints.finite_set.FiniteSet(x_bs_ue, vec),  # x\in{0,1}
                   # x_bs_ue >= 0,
                   x_bs_ue <= bs_status_temp,
                   cp.sum(cp.multiply(rb_req, x_bs_ue), 0, keepdims=True) <= np.ones(shape=[1, const.BS_NUM]),
                   # 原本老的约束条件
                   # cp.sum(x_bs_ue, 1, keepdims=True) <= one_ue,
                   cp.sum(x_bs_ue, 1, keepdims=True) == one_ue,
                   cp.sum(x_bs_ue) >= np.ceil(ue_num * const.UE_REQ),  # 用户速率要求条件
                   ]
    c = const.P_WORK * bs_status + const.P_SLEEP * (const.ONE_BS - bs_status)  # shape =[1,bs_num]
    RB_W = np.tile(const.RB_W, (ue_num, 1))
    f1 = cp.sum(
        cp.multiply(cp.multiply(cp.multiply(rb_req, r_bs_ue), RB_W),
                    x_bs_ue)) * 3600 - ce_temp * const.ETA_GRID * cp.sum(
        cp.maximum((const.RB_PR / const.ETA) * cp.sum(cp.multiply(rb_req, x_bs_ue), 0) + c[0] - p_available[0], 0))
    f2 = function2(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req)
    f2 = f2 + cp.sum(
        cp.multiply(concave_derivation(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req), (x_bs_ue - res_x_bs_ue)))
    # f2 = ce_temp * const.ETA_PV * cp.sum(
    #     cp.minimum((const.RB_PR / const.ETA) * cp.sum(cp.multiply(rb_req, x_bs_ue), 0) + c[0], p_available[0]))
    object = cp.Maximize(f1 - f2)
    problem = cp.Problem(object, constraints)
    problem.solve(solver='GUROBI')
    print('f1-f2的值', problem.value)
    print("optimize_concave的值是：",x_bs_ue.value)
    return x_bs_ue.value


def optimize_concave_none(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req, r_bs_ue):
    bs_num = res_x_bs_ue.shape[1]
    ue_num = res_x_bs_ue.shape[0]
    bs_status_temp = np.tile(bs_status, (ue_num, 1))  # 复制拓展矩阵
    one_ue = np.ones(shape=[ue_num, 1])
    vec = np.array([0, 1])

    x_bs_ue = cp.Variable((ue_num, bs_num))
    constraints = [cp.constraints.finite_set.FiniteSet(x_bs_ue, vec),  # x\in{0,1}
                   # x_bs_ue >= 0,
                   x_bs_ue <= bs_status_temp,
                   cp.sum(cp.multiply(rb_req, x_bs_ue), 0, keepdims=True) <= const.RB_MAX,
                   # cp.sum(x_bs_ue, 1, keepdims=True) <= one_ue,
                   cp.sum(x_bs_ue, 1, keepdims=True) == one_ue,
                   # cp.sum(x_bs_ue) >= np.ceil(ue_num * const.UE_REQ),  # 用户速率要求条件
                   ]
    c = const.P_WORK * bs_status + const.P_SLEEP * (const.ONE_BS - bs_status)  # shape =[1,bs_num]
    RB_W = np.tile(const.RB_W, (ue_num, 1))
    f1 = cp.sum(
        cp.multiply(cp.multiply(cp.multiply(rb_req, r_bs_ue), RB_W),
                    x_bs_ue)) * 3600 - ce_temp * const.ETA_GRID * cp.sum(
        cp.maximum((const.RB_PR / const.ETA) * cp.sum(cp.multiply(rb_req, x_bs_ue), 0) + c[0] - p_available[0], 0))
    f2 = function2(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req)
    f2 = f2 + cp.sum(
        cp.multiply(concave_derivation(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req), (x_bs_ue - res_x_bs_ue)))
    # f2 = ce_temp * const.ETA_PV * cp.sum(
    #     cp.minimum((const.RB_PR / const.ETA) * cp.sum(cp.multiply(rb_req, x_bs_ue), 0) + c[0], p_available[0]))
    object = cp.Maximize(f1 - f2)
    problem = cp.Problem(object, constraints)
    problem.solve(solver='GUROBI')
    print('userf1-f2的值', problem.value)

    return x_bs_ue.value


#  用户关联
def ccp_user(res_x_bs_ue, ce_temp, bs_status, p_available, rb_req, l_bs_ue, h_bs_ue):
    '''
    用户关联：内层CCP算法
    :param res_x_bs_ue:基站用户关联变量
    :param ce_temp:暂时的碳效CE*
    :param bs_status:基站休眠状态
    :param p_available:基站可用能源
    :param l_bs_ue: 基站用户距离
    :param h_bs_ue: 基站用户信道状态
    :return: res_x_bs_ue: 基站用户关联变量
    :return: res_x_bs_ue:基站用户关联变量
    '''
    bs_num = l_bs_ue.shape[1]  # 求取基站数
    ue_num = l_bs_ue.shape[0]

    #  速率R
    one_bs_ue = np.ones(shape=[ue_num, bs_num])
    r_bs_ue = np.log2(one_bs_ue + (const.RB_PR * rb_req * h_bs_ue) / (np.power(l_bs_ue, const.ALPHA) * const.NOISE))

    # bs_status_temp = np.tile(bs_status, (ue_num, 1))  # 复制拓展矩阵

    # 初始无接入
    # res_x_bs_ue = np.zeros([ue_num,bs_num])
    # res_x_bs_ue = nearby_access(bs_pos, ue_pos, bs_status, rb_req)
    converged = False
    x_list = [res_x_bs_ue]
    I = 100
    i = 1
    while not converged and i <= I:
        # print(i)
        res_x_bs_ue_old = res_x_bs_ue
        # 优化凹函数
        res_x_bs_ue = optimize_concave(res_x_bs_ue_old, ce_temp, bs_status, p_available, rb_req, r_bs_ue)
        # 如果为None,说明无法满足需求, 则调用无用户在线率约束的优化方法
        if res_x_bs_ue is None:
            res_x_bs_ue = optimize_concave_none(res_x_bs_ue_old, ce_temp, bs_status, p_available, rb_req, r_bs_ue)
        x_list.append(res_x_bs_ue)

        # 判断是否收敛
        error = np.abs(
            function(res_x_bs_ue_old, ce_temp, bs_status, p_available, rb_req, r_bs_ue) - function(res_x_bs_ue, ce_temp,
                                                                                                   bs_status,
                                                                                                   p_available, rb_req,
                                                                                                   r_bs_ue))
        if error < 1e-05:
            converged = True
            break
        else:
            i = i + 1

    if not converged:
        print('用户关联ccp算法未成功收敛！')
    else:
        print('用户关联ccp算法第' + str(i) + '轮收敛！')
    if function(x_list[-1], ce_temp, bs_status, p_available, rb_req, r_bs_ue) > function(x_list[0], ce_temp, bs_status,
                                                                                         p_available, rb_req,
                                                                                         r_bs_ue) and converged:
        print('用户关联ccp算法有效且成功收敛')
        print('就近接入结果：', f"{function(x_list[0], ce_temp, bs_status, p_available, rb_req, r_bs_ue):.2e}")
        print('用户关联ccp算法结果：', f"{function(x_list[-1], ce_temp, bs_status, p_available, rb_req, r_bs_ue):.2e}")
    else:
        print('用户关联ccp算法无效！')
    return res_x_bs_ue


# np.random.seed(const.SEED)
# bs_num = const.BS_NUM
# ue_num = const.UE_NUM_SEQUENCE[0]  # 10
# bs_pos = uniform2d(bs_num - 1, const.L)
# bs_pos = macro_bs_init(bs_pos, const.L)
# ue_pos = uniform2d(ue_num, const.L)
#
# # 初始化光伏能量
# # p_pv = np.random.poisson(const.PV_A, [1, const.BS_NUM])  # shape=[24, bs_num]
# # s_pv = [const.S_PV_MICRO for i in range(const.BS_NUM - 1)]
# # s_pv.insert(0, const.S_PV_MACRO)
# # p_pv = p_pv * s_pv * const.GAMMA_PV
# p_pv = np.zeros(shape=[1, bs_num])
# # 初始化用户要求和信道状态
# ue_snr_theta = np.random.randint(const.SINR_MIN, const.SINR_MAX, size=[ue_num, 1])
# h_bs_ue = channel_status(bs_num, ue_num)
#
# bs_status = np.array([[1, 1, 1, 1]])
# ce_temp = 232300000000
# l_bs_ue = distance_calculation(bs_pos, ue_pos)  # 测算距离
# _, rb_req = rb_demand(ue_snr_theta, l_bs_ue, h_bs_ue)
# res_x_bs_ue = nearby_access(bs_pos, ue_pos, bs_status, rb_req)
#
# result = ccp_user(res_x_bs_ue, ce_temp, bs_status, p_pv, rb_req, l_bs_ue, h_bs_ue)
#
# print("\n就近接入的CE的最大值",
#       f'{compute_CE(res_x_bs_ue, p_pv, bs_status, h_bs_ue, l_bs_ue, rb_req):.2e}')
# print("\n用户共享优化后CE的最大值",
#       f'{compute_CE(result, p_pv, bs_status, h_bs_ue, l_bs_ue, rb_req):.2e}')
#
# _, _, delta_connected = rate_sum_calculation(result, h_bs_ue, l_bs_ue, rb_req)